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A  MAXIMUM-LIKELIHOOD  APPROACH  TO  IMAGE  SEGMENTATION  BY  TEXTURE 

I.  E.  Bevington  and  R.  M.  Kersereau 


Georgia  Institute  of  Technology 
School  of  Electrical  Engineering 
Atlanta,  Georgia  30312 


ABSTRACT 


This  paper  addresses  the-  problem  of  segment¬ 
ing  an  image  by  texture.  Specifically,  we  are 
concerned  with  estimating  the  trajectory  of  the 
boundary  between  two  regions  characterized  by 
different  two-dimensional  autocorrelation  func¬ 
tions*  An  algorithm  Is  developed  which  is  ah 
approximation  to  a  maximum-likelihood  estimator 
based  on  the  assumption  of  Gaussian  random  fields 
with  known  mean  and  covariance.  Preliminary 
experimental  results  are  given. 


INTRODUCTION 


The  problem  of  segmenting  an  image  Into  more 
or  less  homogeneous  regions  arises  In  the  areas  of 
image  coding  and  Image  understanding.  Segmenta¬ 
tion  may  be  based  on  a  nuaber  of  local  properties, 
one  of  which  Is  texture.  Although  many  defini¬ 
tions  of  texture  are  available,  we  will  for  this 
work  assume  that  a  textured  region  may  be  charac¬ 
terized  by  its  statistical  autocorrelation  proper¬ 
ties. 


The  specific  problem  addressed  In  this  paper 
it  that  of  estimating  the  boundary  between  two 
regions  which  have  different  textures  but  not 
necessarily  different  average  gray  levels.  The 
algorithm  which  Is  presented  assumes  complete 
knowledge  of  the  mean  and  autocorrelation  func¬ 
tions  which  characterise  the  two  textures.  In 
practice,  this  information  would  be  estimated  from 
neighboring  regions  which  are  assumed  homogeneous. 
We  envision  then  that  this  algorithm  might  be  used 
following  an  Initial  coarse  segmentation  proce¬ 
dure,  which  divides  the  image  Into  frames  and 
pre-classlf les  each  frame  as  described  in  |1). 


The  boundary  estimation  algorithm  is  derived 
as  s  maximum-  likelihood  estimator  on  the  assump¬ 
tion  of  Jointly  Caustlan  statistics  for  the  Image 
data.  The  Gaussian  assumption  Is  made  for  mathe¬ 
matical  ‘ractablllty.  Although  actual  Image  data 
may  not  fit  the  Cauaslan  model  very  well,  the 
algorithm  developed  under  this  assumption  has 
parformed  well  in  Initial  testa  with  real  Image 
data.  The  Caustlan  aaaunptlon  leads  to  an  Imple¬ 
mentation  which  involves  two-d loons ional  linear 
prediction.  In  this  sense  our  work  closely  paral¬ 
lels  the  work  of  Quatleri  (2),  in  which  a  statis¬ 
tical  formulation  of  the  problea  of  detecting  an 


anomaly  In  an  otherwise  homogeneous  textured 
region  also  led  to  a  2-D  I.PC  Implementation. 


ESTIMATOR  FORMULATION 


The  2-D  boundary  estimation  problem  la  that 
of  finding  a  contour  C  which  partitions  a  sat  of 
observed  image  samples  (x(m,n))  Into  two  sets, 
where  one  set  is  presumed  to  be  s  realization  of 
the  discrete  rsndoa  field  (Xj(m,n)},  the  other  of 
the  field  {x^Ca.n)}.  We  assume  that  (Xj(m,n))  and 
{x2(m,n)l  are  independent  homogeneous  Gaussian 
random  fields  with  known  mean  and  covariance. 


For  this  paper  we  Impose  the  following  con¬ 
straint  on  the  contour  C  and  the  resulting  parti¬ 
tion  of  the  N  X  N  set  (x(m,n)).  We  require  that  C 
be  characterized  completely  by  the  N-dimenslonal 
vector  L  •  (L(0) ,l( 1 ) , . . . ,l(N-l ) ) ,  such  that 


x< 


(x j(a,n) ,  0<n<L(m) 

* 

x^(s,n),  UaXnOM 


(1) 


where  fK»<N-i . 


That  la,  we  require  that  there  be  only  one 
Hboundary  point"  per  line.  (See,  for  example. 
Figure  1.)  This  requirement  may  at  first  seem 
rather  severe,  but  bear  in  mind  that  we  Intend 
that  the  algorithm  he  applied  only  over  small 
areas  of  an  Image,  and  that  we  have  the  freedom  to 
rotate  the  frame. 


To  formulate  the  maximum  likelihood  estimator 
for  L,  we  need  an  expression  for  the  Joint 
probability  density  function  (pdf)  for  the  set 

(x(m,n))  conditioned  on  L.  This  function  is  most 
easily  expressed  In  vector  notation.  1st  x  be  s 
vector  containing  all  of  the  observed  samples 
(x(iB,n)}.  Now  let  us  partition  x  by  assigning  to 
vector  Xj  those  points  determlned'by  s  particular 
L  to  be  from  <Xj(m,n)),  and  to  x^  those  points 
from  fx  ,(n,n) ) •  The  conditional  pdr  for  x,  which 
we  denote  by  p(x|U  la  the  Joint  density  of  x,  and 
x ,,  conditioned  on  L.  Since  the  random  ?lelda 
Jx  (m,n))  and  (xjm.n))  are  Independent,  thla 
Joint  density  la  simply  the  product  of  the  indi¬ 
vidual  conditional  densities.  That  la, 


P<i| Jr)  “  P(x  ]!l^(jS. 

The  Gaussian  assumption  then  leads  to  (13)) 


(2) 
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Approved  for  public  r-leass: 

distribution unlioi1-'  * 


g jfc  *  -  • 


p(x|L)  •  con«t  •  |  K  j |  “ 1  •  |  K | “ 1 

-  *xJ-7  <5> 

V  *  t'I 


■f  (  x  2  >  T|C  Z  1  (  X  ^ -n  , )]  } 


where 


m,  -  E(x  | 

K!  *  E lilijl  K2  *  EllJip 


Note  that  all  of  the  vectors  and  nitric?*  on  the 
r l ght -hand-side  of  (3)  (even  their  dimensions) 
depend  on  L.  The  explicit  forms  of  these  element* 
also  depend  on  the  order  we  cl»oo*>e  in  assigning 


the  points  of  {x(m,n)}  to  x^  and  x^« 


We  now  consider  a  whitening  transformation  on 
the  observed  data,  which  will  lead  to  an  expres¬ 
sion  for  p(x|L)  which  is  more  easily  evaluated  for 
varying  L.  Let 


«i  -  Vij-V 


(4) 


where  A,  la  such  that 


E<ii*Ii  *  °i 


for  D,  diagonal.  Then  from  (4)  and  (}), 


D1  *  VlAl 


(6) 


Aaauain,  Aj  And  are  Invertible,  this  give. 


Kj1  -  AjDiV, 


(7) 


for  ( x  — ■  .)  .  The  linear  prediction  coefficients 
In  the"Vw©  'matrices  will  depend  on  the  statistics 
of  the  corresponding  processes  and  on  the  ordering 
used  to  form  the  vectors  x  j  anJ  x^ 


Returning  to  p(x|L)  In  (3)  we  take  the 
logarithm,  multiply  by  -2,  and  drop  Che  constant 
term  to  obtain  the  log-1 Ikellhood  function  A ( L) : 


A  (L)  *  1  n  |  K  .  |  4  ln|»2t-  (x  f*  j)’*  '/(*  j-"  j) 

(10) 


(x  P ( x^-«  ^) 


Note  that 


|K;I  -•|A*1|  |d1IIa‘tI 


FI  dj(n) 


on 


(Similarly  for  |k,|.) 

From  (8), (9), (10)',  and  (II), 


A  (L) 


£  fcj(n)/d,(n)  4  ln(dj(a))l 
n  L  J  (12) 

+  £  £ej(n)/d2(n)  4  ln(dj(n))J 


The  aaxlixu.  likelihood  eaclmete  for  L  fa  that 
value  which  minimize.  A (L). 


COMPUTATION  OF  THE  LIKELIHOOD  FUNCTION 


Ualng  (4)  and  (7)  we  can  rewrite  the  firat  half  of 
the  exponent  of  p(x  L)  In  (3)  aa 

(x  j-*,)TK  ( l(Xj-jj)  •  £  e^(n)/dj(n)  (B) 

where  Cj(n)  la  the  nth  element  of  and  d,(n)  ia 
the  nth  diagonal  element  of  D|,  which  la  the  otcan 
aquar.d  value  of  «,(n).  Ualng  almllar  deflnltiona 
for  tj,  Aj,  Oj,  e2(n)  and  d^n)  we  hava 

<*2-»2)1’'2l<i-j'i2)  ‘  £  «  j( n ) /d 2( n  1  (9) 


Now  let  ut  address  the  construction  of  the 
matrices  Aj  and  A^.  All  we  have  required  to  this 
point  la  that  these  operators  whiten  the  vectors 
(Xj*mj)  and  (x and  that  they  be  Invertible. 
One  such  construction  is  aa  full  iws.  (See 
Therrien  14].)  Let  Aj  be  lower  triangular  with  1 's 
on  tha  diagonal.  In  addition,  let  the  non-sero 
off-diagonal  pert  of  each  row  j  be  the  (J-l)st 
order  linear  predictor  for  the  jth  element  of 
(itj-mj)  given  ell  previous  elements.  In  effect, 
we  compute  the  jth  element  of  Cj  by  subtracting 
from  the  Jth  element  of  (xj-mj)  tin*  part  predict¬ 
able  from  ell  elements  k  Tor  which  k<J.  Let  A2 
be  constructed  similarly,  ualng  Iln.nr  predictors 


Wc  now  consider  the  computation  of  the  linear 
prediction  residuals  ej(n)  and  e2(n)  for  a  given 
L.  The  first  step  is  to  establish  an  ordering  for 
the  elements  of  ^  and  x^,  which  we  show  In  Figure 
1.  For  a  given  L,  the  points  assumed  to  be  from 
|xj(n,n))  are  scanned  left  to  right,  top  to 
bottom,  while  those  from  (x  ,(m,n)}  are  a  canoed 
right  to  left,  top  to  bottom.  It  will  be  shown 
that  this  ordering  leads  to  a  simple  llne-by-llne 
procedure  for  estimating  L. 

An  exact  evaluation  of  k  (L)  In  (12)  would 
require  that  the  residuals  «.  and  e2  be  computed 
using  spatially  varying  prediction  masks  for  which 
the  order  increases  with  each  new  point  as  shown 
in  Figure  2.  Also,  the  shape  of  the  masks  for  a 
given  pixel  location  (m,n)  would  depend  on  L(k) 
for  all  k<m  (the  hypothesised  edge-point  locations 
for  all  previous  rows).  In  order  to  arrive  at  an 
implementable  algorithm  we  replace  the  ideal  spa¬ 
tially  varying  masks  with  relatively  small  spa¬ 
tially  invariant  ones  as  shown  in  Figure  2.  There 
are  two  source#  of  error  Introduced  by  this 
approximation.  First,  by  restricting  the  mask 
slta,  we  are  not  in  general  able  to  remove  all  of 
the  correlation  between  the  residual  terms. 
Second,  by  not  allowing  the  shape  of  the  mask  to 
vary  with  Its  location  relative  to  the  hypoth¬ 
esised  boundary,  we  encounter  situations  like  that 
shown  in  Figure  3,  where  part  cf  the  mask  extends 
over  the  boundary,  and  we  In  essence  attempt  to 
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predict  •  unple  from  000  process  using  samples 
assumed  Co  be  part  of  thr  other.  There  are 
obviously  some  tradeoffs  to  be  made  In  selecting 
the  silk  shape  and  at  re. 

The  spatially  Invariant  prediction  mask 
Approx last  ion  leads  to  some  simplification  of 
(12).  because  the  random  fields  (Xj(m,n))  and 
(x2(«(n)}  arc  homogeneous,  the  mean  squared  value 
of  the  prediction  error  is  constant  and  we  can 
replace  dj(n)  and  d-,(n)  by  a  J  and  0*  ,  respec¬ 
tively.  Also,  eaefi  residual  term  e.  or  «2  is  now 
computed  independently  of  L,  *0  the  dependence  of 
A  on  L  is  reduced  to  a  form  which  may  be 
express#?  explicitly  by  the  limits  on  the  summa¬ 
tion*.  Letting  A'(L)  denote  the  approximate 
value  of  A (L)  we  have 

*  N-l-p 

A'U)  -  £  T(m , L(m) )  (13) 

m-p 

where 

k 

T(m,k)  •  £  e^(m.n) /o*  ♦  ln(o^) 
n-p 

N-l-p  ^ 

♦  £  »J(m,n)/05  ♦  ln(o j) 

n»kel 

(The  limits  on  the  aiasmatlons  have  been  chosen 
such  that  we  compute  prediction  residuals  only  to 
within  a  distance  p  of  the  frame  borders,  where  p 
la  the  extent  of  the  predictor  mask.  This 
measure,  taken  to  prevent  the  mask  from  ever 
extending  outside  the  frqme,  slightly  decreases 
the  area  over  which  we  can  search  for  the  con¬ 
tour.) 

As  Implied  In  the  discussion,  the  residuals 

e^m.n)  end  ejCwi")  ere  given  by 

e1(a,n)  •  x(m,n)  -£  £  a,  ,x(is-l  ,n-1) 

1  1 

e2(e,n)  •  *<m,n)  -  £  £  b,  .nfm-i  ,n-J) 

1  1 


Constraints  on  contour  shape  My  be  Imposed 
conveniently  In  the  framework  of  the  expression 
for  A'(L)  in  (13).  In  our  current  implements- 
tlon,  we  first  generate  Tf»,n)  (the  "likelihood 
terrain")  for  all  (m,n)  in  the  search  irea.  The 
next  step  Is  to  search  for  a  path  through  this 
terrain  which  is  in  some  sense  optimal  subject  to 
a  desired  constraint.  We  require  of  course  that 
this  path  touch  exactly  one  point  In  each  line. 
If,  for  example,  the  constraint  la  that  the 
contour  he  a  straight  line,  we  might  systemat¬ 
ically  generate  all  L  which  would  correspond  to  a 
straight  line,  and  evaluate  A  for  each  such  L  by 
summing  the  heights  of  the  terrain  at  points 
touched  hy  the  corresponding  line.  Another  form 
the  constraint  may  take  la  stochastic.  For 
example  we  might  model  the  sequence  of  line-by- 
line  edge  coordinates  as  a  first  order  Markov 
process,  as  is  done  by  Cooper  and  Sung  |5).  In 
this  case  the  optiMl  path  la  found  using  a 
dynamic  programming  procedure. 


PRELIMINARY  RESULTS 

The  results  of  some  early  testing  of  the 
algorithm  appear  in  Figure  4.  We  ahow  results  for 
two  different  32X32  textured  itMges,  one  of  which 
la  repealed  along  the  top  of  the  figure,  the  other 
repeated  along  the  bottom,  the  Images  were  formed 
using  regions  extracted  from  a  multiple-textured 
Image  In  the  USC  data  base.  The  true  edge  (the 
one  used  In  constructing  the  images)  Is  superim¬ 
posed  on  the  first  frame  In  each  row.  The 
difficulty  of  finding  the  edge  by  naked  eye 
(particularly  for  the  second  lmige)  la  Illustrated 
In  the  second  frime  In  each  row,  whure  we  repeat 
the  first  frsme  without  drawing  the  edge.  In  the 
third  frame  we  superimpose  the  boundary  found  by 
our  algorithm  when  we  impose  a  straight-line 
conatraint.  The  fourth  frame  shows  the  result  of 
using  a  crude  stochastic  constraint,  where  the 
edge  coordinate  for  llue  j  la  assumed  to  be  the 
coordinate  on  line  j-1  plus  a  bias  plus  a  uniform 
random  variable  In  the  Interval  (-1,1).  The  bias 
term  is  estimated  hy  first  finding  the  edge 
subject  to  a  straight  line  conaralnt.  The  predic¬ 
tion  masks  used  In  all  instances  were  3X3  quarter 
plane. 


where  (iw)  and  are  the  seta  of  2-0  LPC 
coefficients  for  theJ  processes  (xj(m,n)}  and 
(x Jm.n)),  respectively.  The  LPC  coefficients  are 
dejMndent  on  the  prediction  masks  used  (1). 


EDGE  SHAPE  CONSTRAINTS 

The  likelihood  function  In  (13)  la  minimised 
by  choosing  independently  for  each  row  m  the  value 
L(m)  which  minimises  tvs  correspond! ng  T(m,L(n)). 
It  hot  been  found  in  practice  however  that  much 
better  results  are  obtained  when  we  introduce 
const rslnti  on  the  shape  of  the  boundary  so  that 
coupling  Is  Introduced  into  the  row-by-row  deci¬ 
sions.  The  need  for  shape  constraints  arises 
partially  because  of  errore  Introduced  by  the 
email  spat  tally-invariant  prediction  masks. 


We  have  found  that,  In  general,  the  mure  we 
know  about  the  **d>»r  beforehand,  the  more  accurate 
la  the  resulting  carinate  and  the  more  quickly  it 
ri"  :  found.  Wo  are  currently  studying  Intelll- 
g.  ..  irch  procedut  ea  and  ways  of  specifying  and 
incorporating  a  priori  Information.  We  are  also 
investigating  the  effects  of  prediction  mask  sice 
and  shape  on  the  accuracy  of  the  boundary 
astlMte. 
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